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ABSTRACT 


1.  INTRODUCTION 


The  elements  of  a  theory  for  elastic  composites 
with  a  changing  microstructure,  such  as  microcracking, 
are  reviewed.  This  formulation,  which  uses  internal 
state  variables  and  potentials  like  strain  energy  and 
work,  is  then  illustratee  by  mathematically 
characterizing  and  predicting  the  damage  and 
deformation  behavior  of  particle-filled  rubber  under 
axial  straining  and  confining  pressure.  Next,  a 
micromechanical  model,  which  accounts  for  effects  of 
distributed  microcracks  and  particles  on  overall 
deformation  behavior,  is  described  and  shown  to  be  in 
agreement  with  experimental  findings.  It  is  then 
indicated  how  this  model  may  be  combined  with  the 
potential  theory  to  extend  the  results  from  specimens 
under  axial  straining  and  pressure  to  more  general 
strain  states. 

NOMENCLATURE 


a,b,R 
Aij’  Bij 


'1  j 
E 

fm 

P 


wj 

s* 

uj 
v.  ed 
V.  Vo 
W 

WT,  WjD 

e 


Radii  In  model  composite 
Modu 1 i 

Coefficients  in  dual  energy 

p 

Voung’s  modulus  GPa  (Ib/in*-) 
Thermodynamic  force 
Pressure  GPa  (lb/in^) 
Generalized  displacement 
Generalized  force 
Stiuctural  parameters 
Displacement 
Di latation 
Volume 

Strain  energy 

Total  work  and  dual  work 

Strain 

Poisson's  ratio 
Stress 


In  this  paper  we  use  a  theory  based  on 
thermodynamics  with  internal  state  variables 
(Schapery,  1990a)  to  characterize  and  then  predict 
mechanical  behavior  of  a  highly-filled  rubber  with 
distributed  microcracks.  A  similar  study  of  fiber- 
reinforced  plastic  is  reported  elsewhere  (Schapery, 
1989a).  The  previously  developed  micromechanical 
model  (Schapery,  1986)  serves  to  guide  the 
thermodynamic  formulation  for  general  stress  states 
and  to  identify  the  physical  significance  of  the 
internal  state  variables.  The  effect  of  individual 
particles  and  microcracks  is  smeared  out  in  that  the 
composite  is  represented  as  a  homogeneous  continuum  on 
a  scale  that  is  much  larger  than  particle  and  crack 
sizes.  This  type  of  idealization  is  basic  to  the  many 
studies  in  "continuum  damage  mechanics”;  see,  e.g. 
Krajcinovic  and  Lemaitre  (1987).  However,  the 
underlying  thermodynamic  formulation  (Schapery,  1990a) 
is  nnt  limited  to  this  idealization. 

The  nonlinear  behavior  of  filled  rubber,  such  as 
solid  propellant,  tires  and  numerous  other  commercial 
rubber  products,  has  been  studied  extensively  over 
many  years.  Large  strains  (Swanson,  1983),  vacuole 
(or  crack)  formation  and  growth  (Farris,  1968), 
viscoelasticity  (Schapery,  1982),  and  the  Mullins' 
effect  in  unloading  behavior  (Mullins,  1969)  all 
contribute  to  the  complex  behavior  of  solid 
propellant,  wnich  is  the  material  that  is  used  here  to 
Illustrate  the  damage  theory.  In  this  paper  we 
account  for  only  nonlinearity  of  stress-strain 
behavior  during  loading.  Viscoelasticity  and  the 
Mullins'  effect  may  be  incorporated  using  the  approach 
described  in  other  papers  (Schapery,  1982,  1990a).  In 
contrast  to  plasticity  theory,  nonlinear  behavior 
during  loading  is  not  necessarily  tied  to  residual 
(plastic)  strains  or  other  behavior  associated  with 
unloading;  this  point  is  illustrated  by  Schapery 
(1989a)  for  fibrous  composites. 

In  Section  2  we  summarize  the  primary  results 
from  the  thermodynamic  theory  (Schapery,  1990a)  which 
is  expressed  in  terms  of  an  arbitrary  number  of 
independent  generalized  displacements  and  internal 
state  variables.  Characterization  of  the  behavior  of 
solid  propellant  subjected  to  confining  pressure  and 
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axial  stretching  serves  in  Section  3  to  illustrate  the 
theory  when  one  internal  state  variable  is  sufficient 
to  account  for  the  damage.  An  elementary  model  of 
microcracking  is  used  in  Section  4  to  connect  the 
internal  state  variables  to  the  microstructure.  It 
should  be  noted  that  Sections  3  and  A  describe  and 
expand  upon  earlier  work  by  Schapery  (1986,  1987a). 
In  Section  5  this  model  is  used  to  generalize  results 
for  ax i symmetric  loading  to  a  general  strain  state. 

2.  ELASTIC  BEHAVIOR  WITH  CHANGING  MICROSTRUCTURE 

We  consider  an  elastic  structure  or  material 
whose  thermodynamic  state  is  a  function  of  independent 

generalized  displacements  q  j  ( j  =  1,2 . J)  and 

internal  state  variables  S-,(nr  =  1,2,.  ..M)  as  well  as 
temperature  or  entropy;  inelastic  behavior  arises  from 
changes  in  the  S^.  Generalized  forces  Q.-  are  defined 
in  the  usual  way  in  that 

6W 1  =  Q j 6q j  (j  not  sunned)  (1) 

for  each  virtual  displacement  cg^,  where  aw  is  tne 
virtual  work.  Then,  from  thermodynamics 


Qj-  aw/aq j  (2) 

where  W  is  the  Helmholtz  free  energy  (when  temperature 
is  used  as  an  independent  state  variable)  or  the 
internal  energy  (when  entropy,  instead  of  temperature, 
is  an  independent  variable).  For  brevity,  thermal 
effects  will  not  be  considered  here,  and  therefore  we 
shall  refer  to  W  simply  as  the  strain  energy.  The 
generalized  displacements  qj  may  be,  for  example,  the 
uniform  strains  in  a  material  element  and  Qj  the 
conjugate  stresses,  or  qj  and  Qj  may  be,  respectively, 
the  displacements  and  forces  applied  to  a  structure. 

The  internal  state  variables  serve  to  define 
changes  in  the  structure  such  as  micro-  or  macro¬ 
cracking,  and  are  called  structural  parameters. 
Whenever  any  one  $  »0,  we  specify  as  the  evolution 
law, 


where  Ws  =  Wa(S_)  ’s  a  state  function  of  one  or  more 
Sm;  also,  fm  is  the  thermodynamic  force, 

f=  -  3W/3Sm  (4) 

The  left  side  of  equation  (3)  is  the  available  force 
for  producing  changes  in  Sm,  while  the  right  side  is 
the  required  force.  For  any  specific  set  of  processes 
(i.e.  histories  q.-(t)),  equation  (3)  may  not  be 
satisfied  by  some  Jof  the  parameters;  if  it  is  not, 
those  Sm  will  be  constant.  The  subscript  r  or  p  will 
henceforth  be  used  in  place  of  m  to  designate  the 
parameters  that  change,  which  are  taken  to  be  R  In 
number. 

The  total  work  done  on  the  body  by  Qj  during  an 
actual  process  (i.e.,  a  process  for  whichJ parameters 
change  in  accordance  with  equation  (3)),  starting  at 
some  reference  state,  is  denoted  by  Wj, 

WTs  I  Qjdqj  (5) 

where  the  surmiation  convention  for  repeated  Indices  Is 
used.  From  equations  (2)-(5)  we  find  that  Wj  is  a 
function  of  the  state  (q^,  Sm),  and  is  given  by 

WT=  W  +  W$  (6) 


where  W  =  W.  =  0  in  the  reference  state.  Thus,  Ws  may 
be  interpreted  as  that  portion  of  the  total  worn  Wj 
which  contributes  to  changes  in  the  structure. 

The  second  law  of  thermodynamics  provides  an 
inequality  as  a  constraint  on  the  changes  in  state, 

• 

W$=  T  S'  >  0  (7) 

where  T  is  absolute  temperature  and  S'  is  the  entropy 
production  rate.  Even  if  equation  (3)  is  satisfied 
for  any  one  Sr,  this  inequality  may  not  allow  it  to 
change.  Additionally,  instantaneous  values  of  the  Sr 
are  such  that  they  minimize  the  total  work  .hen  the 
body  passes  through  stable  states;  i.e., 

3WT/aSr=  0  (8) 
(c2W(/3:r3Sp)sSr6Sp  >  0  if  6Sr6Sr  >  0  (9) 

It  should  be  emphasized  that  minimum  total  work 
is  a  necessary  and  sufficient  condition  for 
stability,  In  uontrast,  minimum  W,  rather  than  wT,  is 
equivalent  to  the  stability  condition  for  systems  in 
the  neighborhood  of  equilibrium  states,  according  to 
classical  thermodynamic  theory.  When  §  >  0  and 
equation  (3),  with  f  *  0,  governs  changes  in  S_,  the 
system  is  not  in  tne  neighborhood  of  equilibrium; 
regard’ess  of  how  slowly  processes  occur,  the  total 
entropy  production  Ws/T  between  any  two  states  is  the 
same. 

It  is  observed  that  equation  (3)  represents  R 
equations  for  finding  the  Sr  as  functions  of  q.. 
Then,  Wj  =j  W^(qj  ,Sr (q  j ) , $q )  where  the  S_  are  the 
constant  parameters.  From  equation  (5), 

Qj=  aWq/sq j  (10) 

showing  that  the  body  exhibits  hyperelastic  behavior 
during  the  time  any  particular  set  of  parameters  Sr 
undergoes  change.  Because  the  total  work  is  a 
potential  during  inelastic  (and  elastic)  processes, 
the  incremental  stiffness  matrix  is  symmetric. 
Conversely,  given  that  the  stiffness  matrix  is 
symmetric  when  one  or  more  Sr  change,  then  both 
equations  (3)  and  (10)  follow. 

If  forces  act  on  crack  faces  they  have  to  be 
included  in  the  set  Qj  unless  they  are  due  to  contact 
without  dissipation  (sliding  without  a  shear  traction 
or  no  sliding);  in  the  latter  cases,  the  effect  of 
crack  opening  and  closing  may  be  taken  into  account 
through  the  form  of  the  strain  energy  function. 
Sliding  with  Coulomb  friction,  if  significant,  cannot 
be  accounted  for  through  a  work  potential,  and 
therefore  the  stiffness  matrix  is  not  necessarily 
symmetric  during  processes  involving  crack  face 
sliding.  If,  however,  one  can  use  a  potential  to 
characterize  the  relationship  between  crack-face 
forces  and  relative  displacements  between  crack  faces, 
equation  (10)  may  be  extended  to  this  case  by 
including  this  surface  potential  (which  may  depend  on 
additional  structural  parameters)  in  Wy.  Such  a 
simplification  is  applicable  with  surface  free-energy 
effects  (Schapery,  1990a);  also,  it  was  proposed  by 
Schapery  (1989b)  to  account  for  crack-face  friction  in 
ice  under  compression. 

In  the  next  section  we  illustrate  the  use  of  the 
theory  by  combining  it  with  experimental  results  on 
solid  propellant.  One  structural  parameter  (internal 
state  variable)  S  appears  to  be  sufficient  to 
characterize  the  damage  state  for  the  particular 
axisymmetric  loading  used.  For  this  case  we  may  take 
S  =  Ws  without  loss  in  generality.  Then  from 


2 


(21b) 


(11) 

(12) 

3.  A  SPECIMEN  UNDER  AXIAL  EXTENSION  AND  PRESSURE 

Consider  a  cylindrical  bar  (not  necessarily  with 
a  circular  cross-section)  which  is  subjected  to 
specified  axial  displacement  U  and  all-around  pressure 
p.  As  generalized  displacements  select 

qj-  U/L0  ,  q2=  AV/Vo  (13) 

where  fiV  =  volume  change  (>0  for  an  increase),  L0  » 
initial  axial  length,  and  V0  =  initial  volume. 
Equation  (1)  yields 

Qr  Fi/A0  •  Q2'  -p  (14) 

where  F^  is  the  axial  force  above  that  due  to  the 
pressure;  F,  =  0  wher.  the  only  loading  is  due  to  p. 

Equations  (13)  and  (14)  are  not  limited  to  small 

strains.  Let  us  now  replace  the  generalized  notation 
by  engineering  stress  and  strain  notation,  in  which 

c  5  qr  v  5  q2‘  0  =  Qi  (15' 

Note  that  v  is  the  dilatation. 

With  the  strain  energy  density  W  =  W(<.,v,S), 
equation  (2)  becomes 

o  =  aW/3  e ,  p  =  -aW/av  (16) 

Inasmuch  as  c  and  p  are  given  quantities  in  this 

example,  it  is  helpful  to  rewrite  the  theory  so  that 

they  appear  as  the  independent  variables.  For  this 
let  us  introduce  a  dual  strain  energy  density  W^  = 
Wrf(c  ,p,S)  defined  by 

Wd=  W  +  pv  (17) 

The  total  differential  of  Eg.  (17)  is 


equations  (3)  and  (4),  when  $  *  0, 

3W  _  , 

3S  " 

and  from  equation  (7), 

S  >  0 


aW  .  aW,  3W, 

— ^  dc  +  dp  +  Tr—  dS 
3  £  3  p  35 


=  -jp  dc  +  dv  +  pdv  +  vdp  +  dS 
3c  av  aS 


(18) 


From  equations  (16)  and  (18), 

o  *  awd/3c  ,  v  =  3Wd/ap  ,  aW/ 3$  =  aWd/aS  (19a) 


and 

from 

equat 

ion  (11)  when 

S  *  0, 

aW./aS  =  -1 

d 

(19b) 

The 

total 

dual  work  is 

defined  by 

WTd5  WT+  pv 

(20) 

so 

that 

from  equation  (5) , 

wTd=  Jodc  + 

Jvdp 

(21a) 

and 

then 

using 

equation  (6) 

and  W$  =  S, 

w 


Td 


S 


This  reformulation  of  the  basic  theory,  in  which  the 
independent  variables  are  a  mix  of  strain  and  stress 
variables,  is  covered  in  the  more  general  earlier  work 
by  Schapery  (1990a)  using  an  arbitrary  number  of 
generalized  forces  and  displacements;  it  is  shown  that 
Wjd,  just  as  Wj,  is  a  minimum  with  respect  to  Sm  at 
each  stable  stale. 

Equations  (19)  and  (21)  will  be  used  with  some  of 
the  experimental  data  in  Fig.  1  to  obtain  wd,  and  then 


Fig.  1  Stress  and  dilatation  behavior  of  a  filled 
elastomer  (65  volume  percent)  at  four  con¬ 
fining  pressures.  Experimental  data  from 
Farris  (1968). 


we  will  predict  the  stress  and  dilatation  responses 
not  used  in  the  characterization  process.  In  each 
test  a  constant  pressure  was  applied  and  then  the 
sample  was  strained  at  a  constant  rate  to  failure.  It 
is  to  be  observed  that  the  stress-strain  curves  are 
essentially  independent  of  pressure  when  the  strain  is 
sufficiently  small.  In  the  strain  range  for  which  the 
stress-strain  curves  are  the  same,  the  dilatation  v  is 
essentially  zero.  The  dilatation  in  Fig.  1  is 
actually  that  due  to  the  stress  o  and  does  not  include 
the  very  small  amount  due  to  the  pressure  alone.  At 
the  highest  pressure  there  is  very  little  dilatation, 
and  therefore  we  shall  assume  for  purposes  of  curve 
fitting  that  S-0  (undamaged  material)  in  this  case. 


Determination  of  constitutive  functions 


As  a  start  on  determining  W^,  let  us  expand  it  in 
a  power  series  in  e  and  p  and  drop  terms  which  are 
higher  than  second  order, 

Wd=  Cj«  +  c2P  +  \  cne2+  \  c22p2  +  c12 cp  (22) 

where  the  coefficients  are.  In  general,  functions  of 
S.  From  equation  (19a), 

0  -  C,+  Cne  ♦  C12p,  V  =  C2+  C12c  -  C22p  (23) 

Regardless  of  the  value  of  S,  we  assume  o  =  v  =  0  when 
c  =  p  =  0,  and  therefore  Cj  =  c2  »  0.  Equation  (23) 
does  not  account  for  the  small  amount  of  nonlinearity 
in  the  stress-strain  curve  at  p  =  500.  We  shall 
account  for  this  nonlinearity  by  replacing  c  in 
equation  (22)  by 

f  =  f(c)  e  f  JCo5dc'  (24) 


3 


where  o,  is  the  stress  for  p  =  500  and  E  is  the 
initial  sslope  of  the  stress-strain  curve;  note  tijat 
for  the  linear  case  Et',  and  thus  f  =  e  . 

Equation  (22)  becomes 

W^=  c^f  +  2  c22p  +  ci2eP  (25) 
Equation  (19a)  yields, 

°5 

o  -  c ^ ^  ^  +  c^P  •  v  —  c^t  ^  C22P  (26) 

Inasmuch  as  dilatation  due  to  microcracking  is  never 
negative  and  there  is  no  pressure  dependence  of  the 
stress  when  v=0  (cf.  Fig.  1),  we  set  C22  -  =  0 

whenever  equation  (26)  predicts  v<0. 

The  three  functions  of  S  in  equation  (26),  i.e. 
Cii,  will  be  found  by  using  the  experimental  data  in 
Fig.  1  from  the  tests  at  the  two  lowest  pressures. 
Specifically,  we  shall  fit  v  in  equation  (26)  to  the 
two  dilatation  curves  and  fit  0  to  the  stress-strain 
curve  at  the  lowest  pressure.  However,  this  process 
gives  us  coefficients  that  depend  on  e,  not  S.  In 
order  to  find  their  dependence  on  S,  we  first  solve 
equation  (21)  for  S, 

2  2 

S  =  WTd'  Wd  =  -f  ode  vdp  *  Wd  (27) 

where  state  1  is  c  =  p  =  S  =  0  and  state  2  is  the 
current  strained  and  pressurized  state.  It  is  assumed 
that  equation  (19b)  is  satisfied  when  *  >  0,  and 
therefore  it  gives  us  S  =  S(c,p).  This  means  Wyd  Is 
independent  of  path  when  *  >  0,  so  that  we  may 

integrate  between  state  1  (e=0,  p-T)  and  state  2 

(e,p)  by  going  from  (0,0)  to  (0,p)  and  then  from  (0,p) 
to  (e,p).  The  dilatation  is  (essentially)  zero  when 
c — 0 ,  and  thus  there  is  no  contribution  from  the  vdp 
term.  Equation  (27)  becomes 

S  =  /ode1  -  Wd  (28) 

where  0  is  the  stress  at  pressure  p  and  strain  e1; 
hence,  the  integral  in  equation  (28)  is  simply  the 
area  under  the  stress-strain  curve  in  Fig.  1  for  each 
pressure.  The  dual  energy,  equation  (25),  may  be 

rewritten  as 

Wd=  i(0£  +  vp)  +  ^-(f  -  e)  (29) 

It  is  seen  from  equations  (28)  and  (29)  that  S 

can  be  found  directly  from  the  experimental  data  at 
each  pressure  in  Fig.  1,  except  for  c,j  (which  depends 
on  S).  The  approach  we  used  to  determine  S  was  to 
start  with  cii  =  E  in  equation  (29),  and  then  find.c.y 
in  terms  of  this  S,  say  S'1';  an  improved  S,  say  S'^', 
was,  .then  found  by  using  the  previously  derived  Cji 
(S'1')  1r«  equation  (29).  Further  iterations  produced 
no  further  change  in  S  and  Cjj.  For  each  S,  all  three 
coefficients  Cjj,  Cj2.  c22  were  der’ved  using  equation 
(26)  for  two  pressure  levels.  Namely,  for  pressure  p 
=  15  data  for  both  0  and  v  were  used;  introducing  the 
subscript  a,  we  may  write  for  this  pressure, 

V  cllT  +  c12pa*  V  c12e  +  c22pa  (30) 


Pa  in  a  first  approximation,  c  j  j  and  Cjo  were  obtained 
directly  from  equation  (30)  and  the  data,  and  then 
they  were  expressed  in  terms  of  S  using  equation 
(28).  Finally,  equation  (31)  together  with  S  for  pb  = 
65  yielded  c22.  as  Cy2(S)  was  known  at  this  stage. 
The  coefficients  were  recalculated  by  using  pa  =  15 
and  the  first  approximations  for  C[2  and  C22  in  the 
second  terms  in  equation  (30).  Convergence  was 
achieved  with  three  iterations.  Ihe  iterative  process 
described  earlier  for  S  used  the  third  approximation 
for  c j  ^  with  each  iteration.  It  is  not  really 
necessary  to  use  iterations  to  find  Cj<;  we  could  have 
found  all  three  coefficients  by  simultaneous  solution 
of  equations  (30)  and  (31)  with  pa  =  15,  but  It  would 
have  been  necessary  to  interpolate  with  respect  to  c 
since  the  same  value  of  S  has  to  be  used  in  both 
equations. 

Figure  2  shows  the  coefficients  as  functions  of 
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Fig.  2  Coefficients  in  equation  (26)  for  E  =  400  psi 
(2.76  MPa). 

the  structural  parameter  over  the  S-range  found  from 
the  experimental  data.  Low-order  Chebyschev 

polynomials  were  used  to  obtain  an  analytical 
representation  for  each  coefficient.  Second  order  was 
used  for  Ciy  and  Cj2  and  third  order  was  used  for  c22; 
however,  th  .the  case  of  c,i,  the  expansion  was  in 
terms  of  instead  of  S. 

It  should  be  noted  that  the  results  are  not 
limited  to  small  strains,  in  that  they  are  based  on 
experimental  data  out  to  c  =  0.55.  However,  if  the 
geometrically  linear  expression  for  dilatation  is  used 
to  calculate  transverse  strain  e.,  i.e.  2c.+  e  =  v, 
we  find  the  Poisson's  ratio  for  an  unpressurized  bar 
to  be  simply  related  to  c12, 

v  i  -ct/c  =  (1  -  c12)/2  (32) 

where  v  =  has  been  used;  from  Fig.  2,  v  =  0.5  at 
S  =  0  and  v  =■  0.36  at  S  =  10. 

Prediction  of  the  mechanical  state 

Having  found  the  coefficients  in  W.,  equation 
(25),  we  shall  now  predict,  S,  0,  and  v  by  means  of 
equation  (26)  and  the  growth  equation  for  S,  equation 
(19b),  i.e. 


For  a  second  pressure,  =  65,  use  only  v, 

vb=  C12c  +  c22pb  ^31) 

Recall  that  0-  is  a  function  of  only  e,  and  that  the 
coefficients  are  functions  of  only  S.  8y  neglecting 


1  ^11  f  1  ^22 

2  dS  r  2  dS 


2  .  dC12 
P  +  xT" 


■P  =  -1 


For  each  (c,p)  we  may  solve  equation  (33)  for  S  by  the 
Newton-Raphson  method. 

Figure  1  shows  the  theoretical  predictions  using 
continuous  lines.  Recall  that  the  stress  for  p=15  and 
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dilatation  for  p= 1 5  and  p=65  were  used  in  curve¬ 
fitting  the  coefficients.  The  small  discrepancy 
between  theory  and  experiment  for  these  cases  is  due 
to  the  use  of  low  order  polynomials  for  Cjj;  however, 
use  of  higher  order  polynomials  was  founa  to  produce 
less  smooth  curves  (not  making  a  special  effort  to 
reduce  the  roughness)  and  this,  in  turn,  lead  to 
convergence  problems  in  solving  equation  (33). 
Prediction  of  dilatation  for  p  =  1 65  and  500  is 
excellent,  and  the  stress  predictions  are  quite 
good.  It  should  be  added  that  the  original 
experimental  data  were  not  available  to  us.  Rather  it 
was  necessary  to  take  values  from  the  curves  in 
Farris'  (1968)  publication,  which  is  probably  the  main 
source  of  the  difficulty  in  refining  the  curve  fitting 
cf  c^. 

The  source  of  the  error  in  stress  predictions  for 
p>15  is  at  least  due  in  part  to  the  fact  that  the  data 
do  not  satisfy  exactly  the  uondition  for  existence  of 
a  work  potential,  regardless  of  the  number  of 
structural  parameters  which  are  used.  This  condition 
may  be  found  by  differentiating  the  two  relations  in 
equation  (19a), 


av 

dc 

do 

ap 

(34/ 

Integration  yields. 

.  e 

v  = 

3  r  x  i 
»  o  °d£ 

(35) 

showing  how  dilatation  can  be  predicted  by 
differentiation  of  the  area  under  the  stress-strain 
curves.  This  relationship  is  not  fully  satisfied  by 
the  experimental  data  in  F^g.  1;  in  earlier  work,  the 
existence  of  a  work  potential  was  explicitly  checked 
using  a  less  direct  but  equivalent  relationship  based 
on  equation  (16)  (Farris,  1968  and  Schapery,  1987b). 
The  relative  error  in  equation  (35)  is  about  the  same 
as  in  Fig.  1  for  stress  predictions,  and  thus  the 
latter  error  does  not  appear  to  be  due  to  the  use  of 
only  pr,e  structural  parameter.  Viscoelasticity,  which 
is  neglected  here  and  in  these  earlier  studies,  may  be 
at  least  partly  responsible  for  the  discrepancy,  and 
thus  for  the  difference  between  theory  and  experiment 
in  Fig.  1;  a  method  of  accounting  for  viscoelasticity 
is  discussed  elsewhere  (Schapery,  1990a, b). 

Finally,  it  should  be  noted  that  Farris  (1968) 
employed  p=l 5  to  designate  atmospheric  pressure,  and 
thus  used  absolute  pressure  in  characterizing  the 
material  behavior.  We  examined  the  effect  of  using 
gage  pressure  for  p  (so  that  the  pressures  in  Fig.  1 
are  taken  as  0,50,150,  and  485),  and  found  essentially 
the  same  results  as  shown  for  the  theory  in  Figs.  1 
and  2.  When  the  characterization  based  on  absolute 
pressures  was  used  to  predict  the  stress  for  p  =  0,  at 
e  =0.55  this  stress  was  found  to  be  8  percent  less 
than  the  stress  for  p  =  15. 

4.  A  MODEL  FOR  MICROCRACKING 

In  an  earlier  study  (Schapery,  1986)  developed  a 
mathematical  model  for  axisynmetric  deformation 
behavior  of  a  random  particle-reinforced  rubber  (or 
any  other  relatively  soft,  incompressible  matrix)  with 
microcracking;  changes  in  the  microstructure  were 
assumed  to  be  entirely  due  to  microcracking.  Emphasis 
was  on  the  micromechanics  of  a  material  which  is 
linearly  elastic  for  any  given  state  of  damage, 
although  some  effects  of  intrinsic  material 
nonlinearity,  viscoelasticity  and  microcrack  growth 
were  studied.  Here  we  shall  briefly  describe  the 
microstructural  geometry  and  two  different  methods  of 


using  it  in  predicting  overall  constitutive 
functions.  The  model  will  then  be  used  to  provide  a 
micromechanical  interpretation  of  the  behavior 
discussed  in  Section  3  and  to  guide  the  development  of 
a  procedure  that  accounts  for  damage- induced 
anisotropy  in  more  general  stress  states  than  the 
axisymmetric  state  employed  in  the  characterization 
process. 

For  a  recent  evaluation  of  micromechanical  models 
of  linear  elastic  moduli  without  cracks,  see 
Christensen  (1990).  Included  is  a  study  of  the 
limiting  behavior  of  composites  with  an  incompressible 
matrix  and  a  high  filler  volume  fraction.  (Solid 
propellant  typically  consists  of  65-80*  volume 
fraction  of  essentially  rigid  nancies  in  a  matrix 
with  a  Poisson's  ratio  very  close  to  0.5,  and  thus 
this  limiting  behavior  is  of  interest  to  us.)  For  the 
so-called  generalized  self-consistent  method  (GSCM) 
without  cracks,  Christensen  (1990,  Table  1)  gives  the 
ratio  of  composite-to-matr ’  .  shear  modulus  as  27/16(1- 
v  )  ,  where  v  =  1  is  the  particle  volume  fraction. 
For  consistency  with  Schapery' s  (1986)  result,  the 
numerical  factor  should  be  27/4  instead  of  27/16;  in  a 
communication  with  Dr.  Christensen  it  was  confirmed 
that  27/4  is  the  correct  value. 

Geometry  of  the  microstructure.  Analysis  methods. 

In  one  case  discussed  by  Schapery ( 1986)  the 
random  microstructure  for  a  material  like  solid 
propellant  is  idealized  as  illustrated  in  Fig.  3;  see 


Fig.  3  Geometry  for  the  generalized  self-consistent 
method  with  microcracks. 

Cornwell  and  Schapery  (1975)  and  Schapery  (1986)  for 
SEM  photomicrographs  of  solid  propellant.  Each 
particle  is  assumed  to  be  rigid,  spherical  (with 
radius  a),  and  surrounded  by  a  incompressible  rubber 
shell  of  constant  thickness,  b-a.  A  high  volume 
fraction  of  filler  is  assumed  so  that  (b-aj/b^l. 
Outside  of  this  two-phase  sphere  is  an  effective 
medium  with  the  (unknown)  properties  of  the  random 
particle  composite.  One  or  two  axisymmetric  cracks 
are  in  the  rubber  shell;  their  surfaces  are  concentric 
with  the  particle  surface,  may  be  within  the  rubber  or 
at  the  interface,  and  have  sizes  defined  by  polar 
angles  e  and  e,. 

The0  composite  model  in  Fig.  3  is  subjected  to 
outer  boundary  displacements  at  r>>c  which  correspond 
to  those  for  a  uniform  strain  field.  Namely,  the 
displacement  components  u^  referred  to  an  orthogonal 
set  of  Cartesian  coordinates  X;  are,  for  the  far- 
field, 

ui=  E^.Xj  (i.j  =  1,2,3)  (36) 
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where  c^.=  e..  are  the  components  of  a  spacewise 
constant  '"straWi  tensor;  the  sunmation  convention  is 
used  here  and  in  what  follows  in  that  repeated  indices 
are  to  be  summed  over  their  range.  (This  displacement 
representation  is  valid  for  a  completely  general  state 
of  small  or  large  strains.) 

In  the  generalized  self-consistent  method  (GSCM) 
of  predicting  linear  elastic  moduli,  the  strain  energy 
of  the  body  in  Fig.  3  is  equated  to  that  of  a 
hrmogeneous  body  having  the  unknown  effective 
moduli.  By  doing  this  without  cracks  and  then  with 
cracks,  the  moduli  may  be  found  from  the  set  of 
equations  arising  from  the  strain  energy  equalities. 

In  the  second  method,  the  displacements  in 
equation  (36)  are  imposed  directlv  on  the  surface  r=b, 
and  the  strain  energy  is  calculated.  This  energy  is 
then  equated  to  that  for  a  homogeneous  body  with  the 
unknown  effective  moduli;  the  equality  yields 
expressions  for  finding  the  moduli.  Without  cracks, 
the  second  method  provides  an  upper  bound  to  the  shear 
modulus  for  Hashin's  (1962)  composite  spheres 
assemblage.  He  assumed  that  there  is  a  broad  enough 
distribution  of  radii  b  that  the  composite  spheres  fit 
together  to  form  a  continuum;  we  make  the  same 
assumption  here. 

tliscussion  of  results 

We  have  used  both  methods  for  ax  i  symmetric 
straining,  going  beyond  the  work  reported  by  Schapery 
(1986)  and  obtained  the  same  qualitative  results  and 
similar  quantitative  results  for  the  two  sets  of 
moduli,  all  as  functions  of  the  crack  surfac"  area. 
Figures  4  and  5  give  the  predictions  of  these  two 
methods  for  the  elastic  coefficients  defined  by 


Fig.  4  Elastic  coefficients  from  the  generalized 
self-consistent  method  with  c=b. 
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Fig.  5  Elastic  coefficients  from  the  two-phase  model. 


an  abscissa  that  is  proportional  to  the  growth  of 
surface  area  of  the  cracks  in  Fig.  3, 

S  ;  (cosbQ-  coss)k  (37) 

where  either  one  crack  or  two  equal  cracks  of  angular 
size  e  have  been  used.  It  should  be  recalled  that  the 
structural  parameter  S  used  in  Section  3  is  the  total 
work  input/volume  less  the  strain  energy  density  (cf. 
equation  (6)  with  S=W$).  In  the  micromechanical 
model,  this  work-difference  is  the  work  of  creating 
new  surface  area.  Assuming  the  fracture  work/area  is 
constant,  k  -  G  ,  where  Gc  is  the  (constant)  critical 
energy  release  rate. 

The  angle  e  and  factor  k  were  selected  to 
produce  approximately  the  same  decrease  in  Cjj(S)  as 
shown  in  Fig.  2  for  very  small  S  and  for  S=10, 
respectively;  agreement  between  theory  and  experiment 
for  Cu  at-  small  S  is  best  when  b  is  vanishingly 
small  (-10"  )  for  the  GSCM  and  when  8°  =  0.01  deg  for 
the  two  phase  model,  which  are  the°values  used  for 
Figs.  4  and  5.  By  definition  E  =  Cjj(0).  A  study  of 
the  GSCM  with  c  >  b  showed  that  the  closest  agreement 
with  the  experimental  results  was  achieved  by  using 
c=b  (i.e.  no  uncracked  shell  around  the  inner  two- 
phase  cracked  composite);  only  the  latter  case  is 
shown  in  Fig.  4. 

Consider  now  the  prediction  of  stress  and 
dilatation  by  means  of  equations  (26)  and  (33),  with 
f  =  e  .  The  theory  in  Fig.  6  is  in  quite  good 


Fig.  6  Stress  and  dilatation  prediction  based  on 
the  generalized  self-consistent  method  with 
two  cracks. 


agreement  with  the  experimental  results  of  Fig.  1 
except  for  dilatation  at  intermediate  pressures.  We 
assumed  there  Is  no  gas  pressure  in  the  cracks.  Of 
the  four  cases  in  Fig.  4  and  5,  the  GSCM  with  two 
cracks  agrees  best  with  experimental  data.  The  two- 
phase  model  results.  Fig.  5,  lead  to  the  poorest 
prediction  because  of  the  poor  behavior  of  Cj2  and 
dCjj/dS,  especially  at  small  values  of  S.  For  all 
four  theoretical  cases,  the  magnitude  of  the 
coefficient  c22  Is  well  below  the  experimental  result 
in  Fig.  2.  Consequently,  the  dependence  of  dilatation 
on  pressure  is  not  accurately  predicted.  Referring  to 
equation  (26),  It  Is  seen  that  c„  is  a  plane-strain 
bulk  compliance  in  that  c22  ■  v/p  Tbr  t=0. 

Let  us  now  return  to  equation  (37)  and  the 
results  for  the  two-crack  GSCM  In  order  to  relate  k  to 
Gc  and  estimate  the  value  of  Gc.  The  surface  area  for 
two  equal- length  cracks  at  a  mean  radius  R  is 


A  =  4nR  (1-cos  b) 


(38) 


The  volume  of  the  representative  volume  element  (the 
two-phase  core  In  Fig.  3)  is  approximately  4nRJ/3  If 
(b-a)/b  <<  1,  and  therefore  the  fracture  work/volume 
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IS 


S  =  3(cos>i0-  coss)Gc/R  (39) 

showing  that  k  =  3GC/R.  The  value  of  k  =  226  psi 
provides  the  best  agreement  between  theory  and 
experiment  for  at  S  =  10,  and  was  used  in 
predicting  the  stress  and  dilatation  curves  in  Fig.  6; 
the  predictions  depend  on  dCjj/dS  through  equation 
(33),  which  could  have  lead  to  an  optimum  k  that  is 
different  from  that  for  Cij  alone.  (There  is  also 
weak  dependence  on  the  other  coefficients  and  their 
derivatives  since  p*0).  For  k  -  226  psi  as  well  as 
6---  0  and  a  representative  radius  for  particles  in 
solid  propellant  of  R  =  10  Jin  (25wm),  we  find  Gc  = 
0.08  lb/in  (13  J/nr).  The  matrix  in  the  composite 
studied  here  is  polyurethane  rubber;  it  is  encouraging 
that  the  intrinsic  or  threshold  value  of  Gc  is  0.05 
Ib/in  for  an  unfilled  polyurethane  rubber  studied 
previously  (Schapery,  1975).  Finally,  we  note  that 
6  =  22  deg  at  e  =  0.55  for  the  predictions  In  Fig.  6 
when  p= 1 5;  for  the  case  of  one  crack  a  =  37  deg.  When 
P> 15  the  maximum  crack  angles  are  smaller. 

Although  the  damaged  material  is  transversely 
isotropic,  it  is  of  interest  to  compare  the  ci  .•  to 
those  for  a  fully  isotropic  material.  We  find  1nJthe 
latter  case 


Observe 
impl ies 
moreover,  -Ec 


cir  E’  ci2~  1'2u‘ 

(40) 

-Ec22=  2(2-2v)(l+v)  =  3c12-Cl2 

that  c,2»  0.28  at  S=10  in  Fig. 
for  art^  isotropic  material  -c. 

2,  which 

,E  =  0.76; 

moreover,  -tc22/c,2=  3  at  c12  =  1.  Both  experimental 
(Fig.  2)  and  Theoretical  results  (Figs.  4  and  5)  show 
that  the  plane  strain  bulk  compliance  c?2  is 
considerably  less  than  that  for  an  isotropic 
material.  This  behavior  is  consistent  with  what  one 
expects  physically  because  the  dominant  orientation  of 
microcrack  normals  is  In  the  direction  of  axial 
stretching.  However,  the  microcrack  model,  in  which 
all  cracks  have  an  axis  of  symmetry  In  the  axial 
straying  dirert<cr,  exhibits  greater  ar.is«t«opy 
than  the  experimental  results.  It  Is  likely  that  an 
Improvement  in  the  model  would  be  gained  by  allowing 
for  a  distribution  of  orientations  8  of  the  crack 
symmetry  axis  x,  relative  to  the  global  material  axis 
xi.  Fig.  7.  This  approach  was  outlined  previously 
for  the  two-phase  model  (Schapery,  1987a)  and  could  be 
implemented  using  the  G$CM;  In  the  latter  case,  two- 
phase  cores  with  different  orientations  would  be 


embedded  in  an  effective  medium.  It  would  be 
necessary  to  account  for  locally  nonaxlsymmetric 
loading  on  each  two-phase  core  (cf.  equation  41)). 
Another  generalization  that  should  yield  improved 
values  of  the  plane  strain  bulk  compliance  is  the 
introduction  of  matrix  compressibility.  Whether  or 
not  the  latter  generalization  alone  would  be 
sufficient  to  obtain  realistic  values  of  c??,  at  least 
for  overall  axisymmetric  loading,  is  not  yet  known. 

Crack  initiation  and  early  growth  behavior 

The  predictions  made  so  far  in  this  section  are 
based  on  the  assumption  that  there  are  either  two 
equal  cracks  or  one  crack  in  the  matrix  layer  around 
each  particle  prior  to  loading;  but  they  are  very 
short  (e  <  0.01  deg)  prior  to  loading.  We  have  also 
studied  the  case  of  two,  initially  unequal  cracks. 
With  axial  straining  it  was  found  that  the  energy 
release  rate  of  the  shortest  crack  reached  the 
critical  value  first  and  that  only  it  grows  until  it 
reaches  the  size  of  the  initially  largest  crack; 
without  cracks  the  highest  stresses  are  at  the  top  and 
bottom  of  the  particle  in  Fig.  3.  Then,  with  further 
straining,  both  grow  as  a  pair  of  equal,  stable 
cracks.  Thus,  when  two  very  short  cracks  pre-exist, 
the  only  situation  of  practical  interest  for  most  of 
the  available  straining  range  is  that  of  two  equal- 
size  cracks. 

The  "pre-existing"  cracks  may  be  indeed  that, 
arising  during  processing.  Most,  however,  probably 
develop  from  the  high  triaxial  tensile  stress 
concentrations  between  particles.  Early  experimental 
and  theoretical  work  by  Gent  and  Lindley  (1959)  showed 
that  the  strength  of  unfilled  rubber  in  equal  triaxial 
tension  is  approximately  5E  /6,  where  E  is  the 
Young's  modulus  (at  zero  strain)  of  the  rubber.  This 
value  was  predicted  by  considering  the  growth  of  an 
initial,  arbitrarily  small,  spherical  void  due  to 
remote,  equal  triaxial  tension.  When  the  tension 
reaches  and  then  exceeds  5Er/6  the  void  is  predicted 
to  become  unstable  (due  to  a  geometric  nonlinearity) 
and  grow  without  bound.  As  hypothesized  by  Williams 
and  Schapery  (1965),  when  the  stress  5Er/6  is 
approached  a  Saturn-ring  crack  forms  and  then  grows 
^untii  there  is  sufficient  local  stress  relief  for 
arrest).  The  smaller  the  initial  cavity  is,  the 
closer  the  strength  is  to  5E _/6;  for  cavities  on  the 
order  of  10um  or  less  we  estimate  for  the  material 
studied  above  that  the  strength  Is  close  to  this 
value. 

Recent  work  by  Gent  and  Wang  (1990)  shows  that 
the  triaxial  strength  of  rubber  Is  somewhat  above  this 
vatue  when  the  size  of  Ine  cavity  is  extremely 
small.  The  new  prediction  accounts  better  than 
previously  for  the  high-strain  behavior  of  rubber;  it 
is  believed  there  is  an  error  in  the  energy  analysis, 
stemming  from  the  use  of  Internal  cavity  pressure 
(instead  of  external  tension)  which  acts  on  the 
fracture  surface,  but  the  general  conclusions  appear 
to  be  correct.  Considering  the  very  high  stress  in  a 
real  highly-filled  composite  propellant  with 
irregularly-shaped  particles,  it  is  reasonable  to  view 
the  "pre-existing"  cracks  in  the  above  analysis  as  the 
arrested  Saturn-ring  cracks.  That  the  experimentally 
observed  dilatation  is  less  than  the  prediction  (Fig. 

6)  at  small  strains  may  be,  at  least  in  part,  due  to 
the  lack  of  initial  cracks  until  the  critical  stress 
with  a  value  of  approximately  Er  between  particles  Is 
reached.  r 


Fig.  7  Coordinate  systems. 
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5.  GENERALIZATIONS  EGO  AGuiTRARY  SlKAiN  STATES 

2 

A*  i  symmetric  damage  clP  (Aj2  ~  A11A22^C22 


With  axisymmetric  damage  the  composite  material 
is  transversely  isotropic.  Regardless  of  the 
micromechanica 1  model  employed,  one  may  show  that  the 
strain  energy  density  in  this  case  for  a  linearly 
elastic  material  can  always  be  written  in  the  form 


V  2  |A11V  +  A22ed  +  2A12edv 


44' 


2 

*13 


'23' 


+  A66(v12 


O 


(41) 


where  xj  is  the  axis  of  material  symmetry  and 


V  e22  "  ‘11  ’ 

ed  =  c33*  v/3 

v  =  cn+  t22+  c33  , 

Y 12"  2t 12 

(42) 

y23e  2c23  ’ 

y13;  2 ' 13 

and  the  five  coefficients 
(which  depend  on  the  state 

A.  •  are  elastic  moduli 
of J cracking  or  damage). 

The  strains  e.  .  are  the  tensor  components  introduced 
in  equation  (36j  . 

The  work-conjugate  stresses  are 


c ]2  (Ai2‘  a22/3)c22  (4Z) 


c2Z=  (2 A12/3  -  An-  A22/9)_ 

The  use  of  A^  in  the  micromechanical  model  not 
only  simplifies  theJ analysis,  but  provides  some  moduli 
which  are  practically  independent  of  damage  over  the 
range  studied  in  Section  4.  Namely,  as  6  increased 
from  zero  to  its  maximum  value  at  the  strain  of  0.55, 
we  found  Aj]/E  decreased  smoothly  from  infinity  to 
about  unity,  while  A[2  and  decreased  by  about  25 
and  15  percent,  respectively.  The  remaining  two 
coefficients,  A pp  and  Agg,  changed  by  7  percent  at 
most.  This  result  implies  the  energy  release  rate  for 
each  crack  in  a  two-crack  composite  core, 


ERR 


3  W  r.  3  W 

4.3  0 _  _  __R _ i 

T"D  3(crack  area)  ~  '  3sine  36 


(48) 


is  dominaleo  by  dA../de  unless  the  dilatation  is  very 
small  relative  to  tne  other  strains.  Thus 


ERR 


Rv 


2  dA 


11 


6sine  de 


(49) 
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o 

o 

T 

1 
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v=  3WQ/3V  =  Anv  +  A 

dE  3V3V  A12v  +  A 
s;  3V3V  A66es 

3Wo/3y12=  A6 
13=  3Wo/3t13=  A44t13 
23=  3Wo/3Y?3=  A44y23 


12ed 

22ed 


The  crack  growth  condition  is  ERR  =  Gc;  this  equation 
implies  6  depends  mainly  on  a  scalar  invariant,  the 
dilatation.  We  may  therefore  consider  the  dilatation 
to  he  at  least  a  rough  measure  of  the  magnitude  of 
damage  during  loading;  we  have  verified  this 

(43)  numerically  using  both  methods  in  Section  4  for  the 

case  of  a  specimen  loaded  by  axial  straining  and 
pressure.  Now,  as  A13  scales  wjth  E,  the  initial 
Young's  modulus,  we  may  write  for  v  >  0, 

B  -  function  (REv2/Gc)  (50) 


The  virtual  work  condition  6W'  =  o.6c..,  for  each 
(i.j)  pair,  yields  J  J 

V  3  (°11+  °22+  °33} 

Oj*  033-  (°13+  °22)/2  (44) 

V  (°22-  °n>/2 

where  o.  .  are  the  components  of  the  stress  tensor  in 
the  Xj  coordinate  system.  For  an  isotropic  material 


All=  * 

A44=  A66=  G 


A2z-  3G  * 


ai2=o 


(45) 


where  K  and  G  are  the  bulk  modulus  and  shear  modulus, 
respectively. 


K  =  E/3(l-2v)  ,  G  =  E/2(l+v)  (46) 

In  deriving  moduli  for  the  micromechanical  models 
based  on  the  geometry  in  Fig.  3,  we  found  the  form  of 
W  and  the  associated  strains  In  equations  (41)  and 
(a2)  to  be  especially  convenient.  It  was  the  A^j  In 
these  equations  that  were  developed  first  usingJthe 
two  methods  discussed  in  Section  4.  Then,  by  relating 
the  strains  and  stresses  in  equations  (42)  and  (44)  to 
those  in  equation  (26)  (with  Oc/E  =  c),  we  obtained 
the  c13.  For  reference,  the  reldbionshios  are 


which  is  an  increasing  function  of  the  argument;  when 
v<0,  the  condition  of  constant  damage,  b=0,  should  be 
used  unless  significant  healing  occurs. 

Regardless  of  whether  or. iy  dA../dB  or  the  other 
modulus  derivatives  as  well  are  needed  in  equation 
(48),  the  material  parameters  affect  the  damage 
through  one  dimensionless  ratio,  RE/GC.  A  decrease  In 
radius  R  or  an  Increase  in  Gc  reduces  the  amount  of 
damage. 

The  strain  energy  density  WQ  in  equation  (48)  is 
considered  to  be  that  for  a  composite  with  many 
particles,  each  with  two  cracks  in  the  example  used. 
As  both  W0VQ  and  the  total  crack  surface  area  are 
proportional  to  the  number  of  particles,  there  is  no 
effect  of  this  number  on  ERR.  Similarly,  the 
work/volume  in  equation  (39)  Is  independent  of  the 
number  of  particles.  However,  both  equations  exhibit 
dependence  on  the  radius  R.  With  a  distribution  of 
sizes,  this  radius  has  to  be  interpreted  as  an 
"effective"  value  for  the  simplified  models 
employed.  If  one  wants  to  explicitly  account  for  two 
or  more  different  sizes,  say  Rm(m=l ,2, . . . ) ,  they  have 
to  be  explicitly  Included  by  using  two  or  more  two- 
phase  composite  cores  (particle  plus  rubber  layer)  In 
the  effective  modulus  analysis.  Growth  of  each  of  the 
associated  crack  angles  e  obeys  equation  (48),  with 
ERR  =  G,..  Each  angle  6  ,  or  the  associated  crack  area 
A_,  Is  the  structural  parameter  Sm  of  Section  2  which 
obeys  equation  (3);  if  area  is  used,  f =  (ERR)m- 
However,  according  to  approximate  equation  (50)  for 
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edLfi  i.i,  dl;  angles  depend  pw-nri  ly  on  one  quantity, 
tne  maximum  dilatation  vm.,A,  .his  maximum  value  is  the 
value  of  d i la'  . "  jn 


current 
loading  history.  Th 
structural  parametr 


if  v-0  for  the  entire 
.nere  is  only  one  significant 
we  may  use  as  this  parameter 


tne  maximum  dilatation,  the  total  fracture  work  W$  for 
all  cracks  (*»  in  Section  2),  or  any  other  quantity 
which  is  re.ated  in  a  one-to-one  fashion  to  W$ .  If, 
however,  v  is  used  the  total  work  Wj  will  exhibit 
some  r  ath-dependence  when  I.-O  because  equation  (50)  is 
n  exac*.  Therefore,  it  is  better  to  use  S  =  W$ 
S’  ce  this  choice  guarantees  path- independence  of  work 
-hen  S-o. 


In  some  cases  the  damage  may  be  axisymmetric 
even  if  the  loading  is  not.  Suppose,  for  example, 

<-22  "r,en  l’2'  anrl  s1'  shear  strains 

vanish'.'  Then,  tt  alBne  determines  the  damage 
orientation,  it  will  be  axisymmetric.  In  such  a  case 
the  experimental  results  in  Section  3  can  be  used 
directly  to  derive  Aj  > ,  Ajp,  and  A^,  without  a 
micromechanical  model.  The  remaining  shear  moduli  A44 
and  A^g  have  been  found  from  the  model  to  be 
insensitive  to  damage,  and  thus  we  may  use  for  them 
the  snear  modulus  in  the  undamaged  state.  This 
approach  has  been  used  to  predict  behavior  in  simple 
snear  from  the  results  in  Fig.  2  (Schapery,  1987a). 

In  general,  the  axis  of  isotropy  (x^  in  Fig.  3) 
maj  not  be  parallel  to  an  axis  in  the  coordinate 
system  used  to  represent  the  overall  (applied) 
strains.  In  order  to  allow  for  an  arbitrary 
or ientat icn  use 


Uj  =  ‘ i j  xj  (51> 

where  i  ! are  the  applied  strains  referred  to 
conveniently  defined  coordinates  x I .  The  relationship 
between  the  strains  in  equations  ^42)  and  (51)  is  of 
course  given  by  the  second-order  tenser  transformation 


>J 


m . .  m  .  t . 

)k  jn  kn 


(52) 


where  irij  ■  is  the  cosine  of  trie  angle  between  the  X; 
and  x a«es.  Substitution  of  equation  (52)  into  (42) 
givesJus  k0  in  terms  of  cl^. 

Nnnaxisyimetric  damage 


For  a  general  state  of  overall  stress  and  strain, 
referred  to  the  xl  coordinates, 


’  i  j  "  BijktLki 


(53) 


If  the  damage  is  axisymmetric  and  we  use  xl  =  x,,  the 
associated  moduli  8...  may  be  easily  expVess'ed  in 
terms  of  the  five  A^jJih  equation  (43).  As  81..  is 
a  fourth-order  tensor,  we  may  then  predict  it  any 
orientation  of  xl.  Both  the  GSCM  and  two-phase  model 
discussed  in  Section  4  can  be  used  to  predict  B!., 
with  nonaxisyimietric  damage.  It  may  be  acceptable  to 
assume  each  crack  is  axisymmetric  (thus  simplifying 
the  two-phase  core  analysis)  as  long  as  one  includes 
cores  with  different  symmetry  axis  orientations 
a  ,  *  (cf.  Fig.  7).  If  the  second  method  of  Section 
4  is  used,  which  is  based  on  the  two-phase  composite 
(rb),  the  global  strain  energy  density  will  consist 
of  a  sum  of  strain  energies,  each  having  the  form  of 
equation  (41),  but  with  the  strains  c-.  expressed  in 
terms  of  i j  .  and  the  axis  orientation  tnrough  equation 

^52)  Tho  f  r\r  f  ho  t  utn  _nhaco  r-r\y  d  in  t  ho  if 


The  JW  for  the  two-phase  core  in  the  GSCM  is 


similar  to  that  in  equation  (41),  but  is  expressed  in 
terms  of  local  boundary  displacements  instead  of  the 
overall  strains  (Schapery,  1986). 

If  an  overall  strain  history  is  imposed  for  which 


the  orientation  of  the  maximum  principal  strain  is  not 
fixed,  one  will  have  to  include  composite  core 
orientations  that  reflect  this  history.  Suppose,  for 
example,  the  principal  strain  orientations  do  not 
change,  but  at  one  time  or  another  each  principal 
strain  has  the  largest  value  and  it  is  large  enough  to 
produce  a  family  of  cracks.  In  general,  the  composite 
will  be  orthotropic,  and  the  moduli  will  depend  on  at 
least  three  structural  parameters,  one  a  (or  S)  for 
each  principal  direction.  The  two-phase  core  in  Fig. 
3  is  simple  enough  that  one  can  predict  the  crack  face 
displacements,  and  thus  readily  introduce  constraints 
which  prevent  crack-face  overlap  when  necessary. 

Different  types  of  damage  may  be  introduced.  If 
the  loading  is  like  that  discussed  in  Section  3,  but 
the  axial  strain  is  compressive  instead  of  tensile, 
the  damage  will  be  axisymmetric,  but  the  cra.k  axis 
will  be  normal  to  the  straining  axis.  A  reasonable 
description  for  this  damage  may  be  a  belt-like  crack 
that  encircles  the  particle  about  the  vertical  axis. 
The  strain  energy  density  of  a  composite  core  with 
such  a  crack  can  be  easily  derived  using  the 
previously  developed  perturbation  solution  (Schapery, 
1986). 

CONCLUSIONS 

A  thermodynamically  based  theory  of  mechanical 
state  behavior  with  damage  growth  was  reviewed  and 
then  successfully  applied  to  a  filled  rubber  under 
axisymmetric  loading.  The  analysis  proved  to  be  quite 
simple  because  of  the  path-independence  of  work  (the 
work  potential)  and  the  requirement  of  only  one 
structural  parameter  (or  internal  state  variable). 

Two  different  micromechanical  models  were  then 
used  to  relate  overall  composite  behavior  to 
microcracking  and  associated  material  parameters, 
including  filler  size;  the  amount  of  damage  is  a 
function  of  these  parameters  through  one  dimensionless 
ratio,  RE/G-,  where  R  is  a  characteristic  particle 
size,  E  is  the  initial  Young's  modulus,  and  Gc  is  the 
critical  energy  release  rate  for  microcracks.  The 
model  based  on  the  generalized  self-consistent  method 
(GSCM),  which  is  a  special  case  of  the  thermodynamic 
theory,  agrees  best  with  the  experimental  results. 

In  view  of  these  encouraging  results,  it  is 
believed  desirable  to  extend  the  GSCM  to  predict 
damage  orientation  and  growth  as  well  as  overall 
mechanical  behavior  under  general  strain  states.  With 
damage  localization,  such  as  the  growth  of  a 
macrocrack  In  the  damaged  continuum,  the  theory  should 
still  apply  except  where  there  are  high  strain 
gradients,  such  as  near  the  macrocrack  edge.  Then,  as 
long  as  the  work  potential  theory  Is  applicable  to  the 
continuum,  the  J  integral  will  be  applicable  as  a 
macrocrack  characterizing  parameter  in  many  cases 
(Schapery,  1987b).  Viscoelastic  effects  may  be  easily 
Incorporated  in  the  theory  using  approximations 
described  elsewhere  (Schapery,  1990b). 
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